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Abstract 

In this paper, we present a mathematical and algorithmic framework for the continuation of point clouds 
by persistence diagrams. A key property used in the method is that the persistence map, which assigns 
a persistence diagram to a point cloud, is differentiable. This allows us to apply the Newton-Raphson 
continuation method in this setting. Given an original point cloud P, its persistence diagram D, and a 
target persistence diagram D', we gradually move from D to D', by successively computing intermediate 
point clouds until we finally find a point cloud P' having D' as its persistence diagram. Our method can be 
applied to a wide variety of situations in topological data analysis where it is necessary to solve an inverse 
problem, from persistence diagrams to point cloud data. 
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1. Introduction 


Let P be a finite set of points in given by 

P= {u, I i= l...,M}. (1.1) 

We call P a point cloud, following the convention in 
topological data analysis (TDA) [2 [2] . TDA pro¬ 
vides us tools to study the “shape” of P. Among 
them, persistent homology mi is one of the most 
useful tools, and it is now applied into various prac¬ 
tical applications, e.g., amorphous solids Hi , pro¬ 
teins [7], and sensor networks [8] (see also [1] and 
references therein). 

Persistent homology can be regarded as a collec¬ 
tion of maps, called persistence maps in this paper, 
from P to a finite set for £ = 0,1, •••, in the ex¬ 
tended plane = IRxlR, where K = Ru{oo}. The set 
Di is called persistence diagram and it encodes the 
Gdimensional topological features of P with met¬ 
ric information (precise definitions are given in Sec¬ 


tion 2.3). 
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In many applications, the point cloud P have an 
intricate “shape” or structure, and, in this situa¬ 
tions, persistence diagrams are used to provide the 
“essential” topological features of P. For exam¬ 
ple, in the papers 0 0, the authors study hier¬ 
archical geometric structures in several amorphous 
solids. In such a case, P is given by an atomic 
configuration of an amorphous solid and consists of 
thousands of points in obtained by molecular 
dynamics simulations. It is a difficult task to di¬ 
rectly study the geometry and physical properties 
of the amorphous solid from P due to its immense 
size. Hence, a key observation of their work is that 
the persistence diagrams of the atomic configura¬ 
tions can capture essential geometric information of 
the amorphous solids. From this significant prop¬ 
erty, using persistence diagrams they obtain various 
physical properties of the solid, such as, finding new 
hierarchical ring structures, decompositions of first 
sharp diffraction peaks, mechanical responses, etc. 

Figure shows a schematic representation of 
Di for silica glass, P, studied in 0 (this corre¬ 
sponds to Figure 1 in that paper). They show 
that the presence of curves in Di precisely distin¬ 
guishes the amorphous state from liquid and crys¬ 
talline states. It means that the normal direc- 
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tions to these curves express characteristic geomet¬ 
ric constraints on atomic configurations of amor¬ 
phous states. Therefore, changing Di => D'^ along 
a normal direction (e.g., black to red in Figure 
and tracing the corresponding deformation P ^ P' 
in the atomic configurations clarify the geometric 
origin of rigidity in the amorphous solid, which is 
currently an important problem in physics and ma¬ 
terial sciences. For this purpose, we need to solve 
an inverse problem: given D'^ find P' in some ap¬ 
propriate setting, and this is the main subject of 
this paper. 



Figure 1: Schematic representation of persistence diagrams. 
Black: Di for silica glass (see for details). Red: A target 
persistence diagram D'^ to study the geometric constraints 
generating the curve Cp. 

In this paper, we present a mathematical and al¬ 
gorithmic framework for solving inverse problems 
of persistence diagrams. Our method is based on 
the continuation method laiiQ], which was origi¬ 
nally developed in numerical bifurcation theory of 
dynamical systems, applied to the setting of per¬ 
sistent homology. More precisely, given a known 
correspondence between a point cloud P and its 
persistence diagram Di, and a target persistence 
diagram D'^, we develop a method to obtain a 
point cloud P' which have D'f^ as its persistence 
diagram. We first represent the persistence dia¬ 
grams as points in some Euclidean space and di¬ 
vide the line segment into small segments 

= D',. Then, 

we solve an implicit equation defined by the per¬ 
sistence diagram for each small segment, us¬ 
ing the Newton-Raphson method, to obtain a new 
point cloud having as its persistence dia¬ 
gram. By successively applying this procedure, we 
finally obtain a desired point cloud P' = . 

We remark that the inverse problem from a per¬ 
sistence diagram Di to a point cloud P is not 


well-posed in general. Namely, it is possible to 
have multiple point clouds giving the same 
(non-uniqueness). Furthermore, the target persis¬ 
tence diagram may not be located in the image 
of the persistence map (non-existence). Our ap¬ 
proach to these issues is as follows: Regarding non¬ 
uniqueness, at each step of the continuation we try 
to find the point cloud closest to the point cloud on 
the previous step of the continuation. This assigns a 
minimality condition on the Euclidean norm of the 
difference for persistence diagrams, and provides 
the uniqueness property. Furthermore, this mini¬ 
mality condition is reasonable in the practical ap¬ 
plications mentioned above, since our input atomic 
configuration is usually realized as a minimum point 
of a certain energy landscape, and hence, finding 
the closest atomic configuration to the minimum 
point is a natural choice. Regarding non-existence, 
we take a practical approach. Namely, we try to 
apply our continuation method and, if the com¬ 
putation is successful we conclude that our target 
persistence diagram is in the image of the persis¬ 
tence map. If not, we investigate the reason for 
non-convergence of the Newton-Raphson method, 
which could be due to non-existence, the presence of 
zero singular values, etc (see Section 4.5). We note 
that it is a very challenging mathematical problem 
to study the image of the persistence map. 

This paper is organized as follows. The funda¬ 
mental concepts, such as simplicial complex mod¬ 
els used to represent the point clouds and persistent 
homology, are introduced in SectionIn Section]^ 
local properties of persistence maps, especially dif¬ 
ferentiability, are studied in detail. Section]^ is the 
core of the paper and is devoted to developing the 
continuation method of point clouds using persis¬ 
tence diagrams. In Section we show some com¬ 
putational examples of the proposed method. Fi¬ 
nally, in Section]^ we conclude with a list of future 
improvements to our continuation method. 


2. Simplicial Complex Model and Persistent 
Homology 

2.1. Simplicial Complex Models 

Let V be a finite set, and E be a collection of 
subsets of V. A simplicial complex is defined by a 
pair (V, S) satisfying (i) {v} 6 E for all u € V and 
(ii) tr 6 E and t c a imply t € E. An element cr e E 
with \a\ = ^ + 1 is called an (.-simplex. 
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Let P be a point cloud (1.11 in For 0 < r < oo, 


we refer to the open ball with radius r as an r-ball, 
and denote it, with center Ui, by 


Br € 


\\x - Mill < r}. 


The Vietoris-Rips complex VR(P, r) of P with 
radius r is defined as the simplicial complex (P, S) 
where the set S of simplices is determined by 


Cr € E Pr,(Ms) n P^(Mt) =#= 0, VMs,M(6(T. 


The definition of the Vietoris-Rips complex de¬ 
pends only on the distances of all pairs in P. Hence 
the Vietoris-Rips complex has an advantage that it 
is computable even if L is large. 

The alpha complex Alp(P, r) pi ITT] is another 
simplicial complex model of P defined by using the 
set of r-balls Br{ui). A significant property of the 
alpha complex is the homotopy equivalence 

M 

\jBr{ui) ^ |Alp(P,r)|, 

i=l 

where |Alp(P, r)| is a geometric realization of 
Alp(P, r). Because of this property, the alpha com¬ 
plex is widely used in practical applications to an¬ 
alyze topological features in P. We note that the 
Vietoris-Rips complex does not satisfy this property 
in general. 

Fast software for computing alpha complexes in 
dimensions P = 2,3 is available, e.g., m- In this 
paper, we use alpha complexes for L = 3, while the 
case L = 2 can be similarly treated. 

For 0 < a < oo, an a-ball B is called P-empty 
if P n P = 0 . For £ = 0,1,2,3, let Deb(P) be the 
set of Psimplices a c P such that there exists an 
P-empty open ball B with dB n P = a. The De¬ 
launay triangulation Del(P) of P is the simplicial 
complex whose simplices are given by Deb(P) for 
£ = 0,1,2,3. 

The three dimensional alpha complex is dehned 
as a subcomplex of the Delaunay triangulation 
Del(P). For each a € Deb(P), let B„ be the small¬ 
est open ball with a c dBu, and pa be the radius 
of Per- Let us dehne Gq.q = P, and Gg^a to be 
the set of Psimplices a e Deb(P) such that Po- 
is P-empty and < a. A simplex in Ua 
is called an attaching Psimplex. The alpha com¬ 
plex Alp(P, a) is defined as the simplicial com¬ 
plex whose simplices are given by Gg^a and their 
faces for i = 0,1,2,3. From this definition, the 
alpha complex Alp(P, a) is a subcomplex of the 


Delaunay triangulation Del(P), and we have that 
Del(P) = Ua Alp(P, a). 

We note that both simplicial complex models de¬ 
fine filtrations of finite type 

VR(P) = (VR(P,r)),,R,„, 

AlpC-P) = (Alp(P,r))^sR^o, 

where M>o is the set of nonnegative reals. Namely, 
we have that VR(P, r) c VR(P, s) and Alp(P, r) c 
Alp(P, s) for r < s and VR(P, r) = VR(P, P) and 
Alp(P, r) = Alp(P, P) for r > P, for a sufficiently 
large P (called saturation time). The radius pa¬ 
rameter r is also called time in this paper, following 
the convention used in persistent homology. 

2.2. General Position 

Let us treat a point cloud P as an ordered set 
induced by the index i = 1,..., M. Then, we can 
assign a single variable u = (mi, ..., um) 6 K” to P, 
where n = LM. Conversely, from a point u 6 K", an 
ordered subset P in with |P| = M can be con¬ 
structed. We explicitly denote this correspondence 
by m(P) and P{u), if necessary, and identify them 
in the following. 

Let X = (V’’)j.gR^j, be a Vietoris-Rips or an alpha 
filtration. For each simplex a 6 V^, we can assign 
its birth radius r^. in the filtration X by the infimum 
radius r satisfying a € V’’. 

In the Vietoris-Rips filtration VR(P), the birth 
radius r^ of a simplex tr = {m^o ,..., m^^ } is a function 
of Ua = {ui„ ,..., Mi J given by 

1 

Xu = - max Mi. - Mi . 

2 0<a<b<g'' " 

We call an edge {ui^,Ui^} that attains the above 
maximum an attaching edge of a. 

Definition 2.1. A configuration u € K" is said to 
be in Vietoris-Rips general position if the following 
conditions are satisfied: 

(i) Mi + Uj for any i + j, 

(ii) Xr + Xr' for any attaching edges t + t'. 

The open set consisting of the points in Vietoris- 
Rips general position is denoted by PyR- 

In the alpha hltration Alp(P), each simplex a 
appears as either an attaching simplex or a simplex 
attached by some attaching simplex t o a. In the 
latter case, the birth radius r„ is given by = r^. 
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Definition 2.2. A configuration u e K" is said to 
be in alpha general position if the following condi¬ 
tions are satisfied: 

(i) u is in general position in the sense of m- 

(ii) rr + Tt' for any attaching simplices t + t' ex¬ 
cept 0-simplices. 


For a simplex cr, let us dehne 





(cr), r = r„, 
0 , r + rc. 


We note that = {((cr)) | a € X^} forms a basis of 
Ci{X). The boundary map di : Ci{X) -> Ct-i{X) 
is defined by linear extension of 


The open set consisting of the points in alpha gen¬ 
eral position is denoted by LAip- 

We note that, in both Vietoris-Rips and alpha fil- 
trations, the condition (ii) implies that an attaching 
simplex is uniquely determined by its birth radius. 

2.3. Persistent Homology 

We briefly review the definition of persistent ho¬ 
mology as a graded module on a monoid ring. Let 
X = (A’')r6R>o be a filtration of finite type with a 
saturation time R. For each X'', let us denote by 
X^ the set of Asimplices in W'". In the following, we 
fix an orientation for each simplex a = {uig ,..., } 

hy if) <■■■< ii, and denote the oriented simplex by 
(cr) = {uig...Ui,). 

Let fc be a field, and let us treat IR>o with a 
monoid structure induced by the addition +. Let 
fc[IR>o] be a monoid ring. That is, fc[]R>o] is a vec¬ 
tor space of formal linear combinations of elements 
of M>o equipped with a ring structure 

(airi)• (a2r2) = {aia2){ri + r2) 

for 01,02 ^ k and ri, r 2 € K>o. 

In the following, the elements in fc[IR>o] are ex¬ 
pressed by linear combinations of formal monomials 
az^, where a e fc, r € K>o, and z is an indeterminate. 
Then, the product of two elements are defined by 
linear extension of 


„ ..ri ^ „ .ri+ra 

QjIZ • (I2Z — Cl\(X2Z 

Let us denote by Ce{X'^) the fc-vector space 
spanned by the Asimplices in A^. The £-th chain 
group Ci{X) is defined as a graded module on the 
monoid ring fc[M>o] by taking a direct sum 

CeiX) = 0 CeiX^) = {(c,) | c, e 

where the action of a monomial z'* on Ci{X) is 
given by the right shift operator 

z'* • (cr) = (cj.) with cj, = 


I Cr-„ r>s, 
I 0, r < s. 


3=0 

where {oj) = {uig... uTj ■ ■. Uif.) is a face of (cr) = 
{uig.. .Uig), and o' means the removal of the vertex 

a. 

The cycle group Zg{X) and the boundary group 
Bi(X) \n Cf,{X) are defined by 

Z^X) = kerdt, Be{X) = imcff+i. 

It follows from dc o di+i = 0 that we have Bg{X) c 
Zi{X). Then, the Ath persistent homology is de¬ 
fined by 

HdX) = ZdX)/BdX). 

The following theorem is known as the structure 
theorem of persistent homology. 

Theorem 2.3 ([4]). There uniquely exist indices 
s,t & Z>o and {bi ,di) 6 for i = 1,..., s, with 

bi < di and bi € IR>o, for i = s + 1,..., s + t, .such that 
the following isomorphism holds 

S / / \ S + t 

HdX)^®Uz^') ©(z^-), (2.1) 

where (z“) expresses an ideal in fc[IR>o] generated 
by the monomial z“. When s or t is zero, the cor¬ 
responding direct sum is ignored. 

The Ath persistence diagram Dg{X) of X is de¬ 
fined as the multiset in determined from the 
decomposition by 

De{X) = {{b„dy\t=l,...,s + t}, (2.2) 

where di = +00 for i = s+1,..., s+t. The pair (bi, di) 
is called a birth-death pair in the Ath persistence 
diagram, and bi, di are called, respectively, the birth 
and death times of the pair. 


3. Local Properties of the Persistence Map 

3.1. The Persistence Map 


Let Di{Xp) be the persistence diagram (2.2) of 
a hltration Xp constructed from a finite set P. By 
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choosing the birth and death times in (2.2) that are 
finite, we can express D^{Xp) as a point 

V — (, d \,..., 6s, dg , 6s+ i,'--;6s+£)cM , 


where m = 2s + t. Then, recalling the identifica¬ 
tion of P and u 6 K", we can regard the persistent 
homology as giving a single correspondence 

W^3u<-^veW^. (3.1) 


In this section, we define an appropriate open set 
O c K" such that this single correspondence is ex¬ 
tended to a map /: O c K" M’" computing per¬ 
sistence diagrams with f{u) = v. 

It should be noted that the dimension m may 
change for a different choice of u € M”. For extend¬ 
ing the single correspondence to a map into M™, we 
use a recent result in [13]. Let us first recall some 
definitions. 

For a metric space {X,dx), the Hausdorjf dis¬ 
tance dn between two subsets A,BcX is defined 

by 


d}i{A,B) = max{sup(ix(a,il),sup(ijc(6. A)}, 

aeA b^B 


Lemma 3.1. Let D be a persistence diagram and 
6 = do^{D,A). If db{D,D') < S, then \D\ < \D'\. 
Furthermore, if d\,{D,D') < e < 5/2, then \D\ = 
\Tg{D')\. 

Proof. Suppose |Z?| > |ZI'|. Then, there exists a 
point in D which is mapped to the diagonal A 
for any bijection 7 : II ^ D'. This leads to 

5 < dh{D,D'), implying the first statement. 

For the proof of the second statement, let 7 be a 
bijection j D ^ D' such that 

dh{D,D') < sup||p- 7 (p)||oo < e. 

p^D 

Suppose dcx,{j{D), A) < e. Then, there exists pe D 
such that doo{j{p),A) < e. On the other hand, we 
have doo{p,j{p)) < sup^,^ ||p - 7 (p)||oo < £■ This 
leads to the contradiction 

6 < doo (p, A) < doc (p, 7 (p)) + dco ( 7 (p), A) < 2e < 5. 

Hence, we have doo{"f{D),A) > e. Moreover, we 
have Tg{D' \ 7 ( 11 )) = 0, otherwise it gives the con¬ 
tradiction ||p-7(p)||oo ^ £• These two prop¬ 
erties show that |ZI| = \Tg(D')\. □ 


where dx{a,B) = inf^eB dx(a, 6) and dx{h,A) is 
defined symmetrically. The Gromov-Hausdorff dis¬ 
tance dcH between two metric spaces {X,dx) and 
{Y,dY) is defined by 


Now let us apply the result in [T3]. Let P and 
P' be two point clouds in with |P| = \P'\ = M. 
Then, since P and P’ are totally bounded, the in¬ 
equalities 


dGii{X,Y) = inf dH(/x^z(Ar),pF^z(F)), 
f,g 

where fx^z and gy^z denote isometric embed¬ 
dings of X and Y into a metric space {Z,dz), 
respectively, and the Hausdorff distance between 
fx^z{X) and gy^ziY) is measured using the met¬ 
ric dz- 

We also recall that the bottleneck distance db be¬ 
tween two persistence diagrams D and D' is defined 
by 

dh{D,D') =infsup||p- 7 (p)||oo, 

where II is a multiset consisting of the points in D 
and the points on the diagonal A with multiplicity 
+ 00 , and 7 is a bijection between D and D'. Here, 
we define the norm ||p||oo = d<x,(p, 0) by the distance 
dco{p,q) = max{|pi - qi\, |p 2 - 92 !} on 
For a multiset II c let 


db{Di{Xp),Di{Xp,)) < dGH{P,P') (3.2) 

for Vietoris-Rips filtrations and 

d^{D,{Xp),Di{Xp,)) <dii{P,P') (3.3) 

for alpha filtrations hold by Theorem 5.2 and The¬ 
orem 5.6 in m, respectively. 

Let 5 = doo{Di{Xp),A) and, for e < 6/2, let us 
set 


= {u' 6 K" I dGH{P{u),P\u')) < e} 
for Vietoris-Rips filtrations, and 

0: = {u'€R^\d^{P(u),P'{u'))<e} 


for alpha filtrations. Then, it follows from 
Lemmathat \Di{Xp(^u))\ = \Te{Di{Xp,(^u')))\ for 
any u' e O’/. Hence, given a single persistence cor¬ 
respondence {u,v), we can define a map 


Tg{D) = {p € II I dco(p. A) > e} 


$ : O: 9 u' ^ T,(iI,(A’B7„o)) 6 K™. 


be the multiset defined by an e-truncation of D from At the end of this section, we show that O)) is an 

the diagonal. open set. We first prove the following lemma. 
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Lemma 3.2. Let P = {xi, ... ,xm} and Q = 
{Ui, ■ ■ ■ tVm} be two point clouds in Set u = 

u{P) e K" and w = w{Q) € K". Then, 

dGii{P,Q) < dii{P,Q) < dRr^{u,w). 

Proof. By definition of the Gromov-Hausdorff dis¬ 
tance, we have 

<^Gh(-P, Q) < dii{P, Q) 

= max{max(iRi {xi, Q); {P^ Hi)}- 


where g constructs the simplicial complex filtration, 
and h computes the persistence diagram. We ex¬ 
tend the decomposition f = ho g of this single cor¬ 
respondence to that of a map /. To this aim, we 
need to construct a proper subset in in which 
the set W is invariant. 

In the case of Vietoris-Rips filtration, W is given 
by all the faces of the (M- l)-simplex for |P| = M, 
independently of the configuration u. Thus, the 
decomposition / = h o g of the correspondence is 
extended naturally to the map 


Without loss of generality, we assume that driiP, Q) 
is achieved by Xi. Then, 

dB.{P,Q) = mindRi,(a:i,%) < ds_L{x^,yf) 

3 


- ,= dM^{u,w), 


and this completes the proof. □ 

Proposition 3.3. is an open set in K". 

Proof. We consider the case of Rips filtrations. Let 
us choose u' e and set P' = P'{u') c . We 
also take a positive real number r > 0 such that 

dGii{P,P') + T<e. 


We claim that the r-open neighborhood of u' 
in K" is contained in For any u” 6 K" 

with dRn(u',u") < T, we have dGii{P',P”) ^ 
dRr.(u' 


') < T from Lemma 3.2 Then, it follows 


from the triangle inequality that 


dGH(P, P") < dGH(P, P') + dGH(P', P") 

< dGH{P,P') + T<e. 


This proves the claim and, since u' is arbitrary, this 
concludes that is an open set in K". We can 
prove the case for alpha filtrations by replacing cZgh 
by dn- 

□ 


3.2. Decomposition of Persistence Map 

Let X = (X’’)j.6Kso be a Vietoris-Rips or alpha 
filtration with a saturation time R, and let us set 
Wi = and W = The correspondence 

(3.1) is decomposed into two parts: 


OvR ^ A K™, 

where OvR = n Uvr, g ■ Ovr 9 u i —r = 
{ra)crsw € and h ■ R‘^ b r i—*■ v e K™ with 
d=\W\. 

On the other hand, in the case of alpha filtra¬ 
tion, note that the set W can generally change de¬ 
pending on the configuration u. Recall that the 
alpha complex Alp(P(it),a) is a subcomplex of 
the Delaunay complex Del(P(u)) and Del(P(M)) = 
Uct Alp(P(M),a). Hence, the set W is given by 
Del(P(M)) in this case. 

For a configuration u satisfying the general posi¬ 
tion assumption, the Delaunay complex Del(P(M)) 
is stable with respect to small perturbations. 
Namely, we can construct an open neighborhood 
Ou c R" of u such that the Delaunay complex is 
invariant in this neighborhood, i.e., Del(P(M)) = 
Del(P'(M')) for all u' 6 0„. Therefore, by setting 
Oaip = n Ou n Paip, we can extend the decom¬ 
position of the map / as 

Oaip ^ K™. 

In [13], the authors study explicit bounds on the 
perturbations of u for which the Delaunay complex 
Del(P) is invariant. 

Remark 3.4. Precisely speaking, h is defined on the 
image of g. 

Remark 3.5. When we consider £-th persistence di¬ 
agram D^(T’), it is sufficient to deal only with H4_i, 
Wi, and Wi+i. 


3.3. Smoothness of Persistence Map 

3.3.1. Vietoris-Rips filtration g 

Lemma 3.6. On the open set OvRj the map g is 

of class C°°. 


r = ira)aiW 


V, 
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Proof. For a simplex a = {ui^,... ,Ui^}, let 

{ui^,Ui^} be the attaching edge, i.e., 



From the assumption of the Vietoris-Rips general 
position, is continuously differentiable on OvR 
whose entries in the lo-th and the i{,-th coordinates 
are given by 

dr a 1 - Mil, dr^r 1 Mil, - '^ig 

2 dui^ 2 Hwi,, ^^all 

and zero otherwise. It is obvious from the same 
argument that r^ is continuously differentiable ar¬ 
bitrarily many times. □ 

Remark 3.7. It follows from the proof of Lemma [3.6| 
that breaking the Vietoris-Rips general position as¬ 
sumption immediately makes that map g loose its 
differentiability. 

3.3.2. Alpha filtration g 

Lemma 3.8. On the open set Oaip, the map g is 
of class C°°. 

Proof. Let cr be a simplex in the alpha filtration 
Alp(P) and r be its attaching simplex. Then, it 
follows from the definition of Alp(P) that the birth 
radius r^ is given by the radius pr of the smallest 
circumsphere of r. 

Let us denote by Ui = (Mi,i,Mi, 2 jthe 
coordinate of each point Ui. Let Ui^ and Ui ^4 be 
given by Ui^ = I and respectively. 

We define the determinant 

Uiiji ■■■ 

Wi2j2 ■■■ '^i2,3k 

'^ik,h '^ik,32 ■■■ '^ik,3k . 

Then, the formulas for for |r| = 2,3,4 are given 
in [TT] as 

T = {ui,Uj} : 

2 + {Mgy 

Pr ^ 

T = {ui,Uj,Uk} ■■ 

T = {ui,Uj,Uk,ui} : 

2 + + 4 • Mf2gg • Mlili 


3l32...3k 


It follows from the definition of the alpha general 
position and the above formula that g is of class 
on Oaip- □ 


Remark 3.9. It follows from the proof of Lemma 3.8 


that breaking the alpha general position assump¬ 
tion immediately makes that map g loose its differ¬ 
entiability. 


3.3.3. The map h 

Let us next study the map h and its differentia¬ 
bility. The map h can be computed by Algorithm 
below, which consists of three parts and is based on 
m and [3]. Each procedure is presented in what 
follows. 


Algorithm 1 Compute Persistence Data from W 
and d _ 

procedure ComputePersistenceData(IF, d) 
B ^ BoundaryMatrix(VF, d) 

P ^ PERSISTENCELEFTRlGHT(i3) 
return PersistenceData(P) 


BoundaryMatrix. For a matrix A = (Ay), let us 
denote its j-th column of A by Aj , and for a non¬ 
zero column Aj , set pivot(Aj) := max{i | Ay + 0}, 
called the pivot index. 

Let us order the set W of all simplices, (Ti < 
(72 < ••• < ctk, by the lexicographical order of 
(Ter,dimer) € IR>o x Z. If two (or more) simplices 
(7, a' appear at the same birth radius with the same 
dimension, we order them by an appropriate rule. 

In this order, a subset {ai,---,ak} for any fc is a 
subcomplex of W in both the Vietoris-Rips and the 
alpha filt rat ions. 

Let C be the vector space spanned by 
{ai,... ,afc}. A matrix representation B = (Bij) 
of the boundary map 9 : C ->■ C is constructed in 
such a way that for an Asimplex ai = {uig,... ,Ui^} 
with io < ■ ■ ■ < ii, the (i,j)-entry is given by 

Jjj — i (~1) J ~ {^*0 I ■ • ■ J ^ik I ■ • • J } J 

I 0, otherwise. 

PersistenceLeftRight. a column operation of 
the form Aj <- Aj + AA^ is called a left-to-right 
operation if A: < j. We call a matrix A' derived 
from A if A' can be transformed from A by left-to- 
right operations. We call a matrix A' reduced if no 
two non-zero columns have the same pivot index. 
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If a reduced matrix A' is derived from A, we call it 
a reduction of A. In this case, we define 

Pa' = \A'j + Q and i = pivot(^')}, 

Ea' = {i \ I < i < KA ^ i' and i + j' for V(i', j') 6 Pa'} 

For a reduction B' of the matrix representa¬ 
tion B, it follows from [4] that we can find a ba¬ 
sis {wi,... ,wk} of C satisfying: (i) the subspace 
spanned by {cti, ... ,(Tfe} is equal to the subspace 
spanned by {uii ,... ,Wk} for any k, and (ii) dwj is 
given by 

By, = I (bi)e-PB' 

^ I 0, otherwise 

Algorithm (a modification of Algorithm 1 in 
m) shows a simple algorithm to compute Pb' of 
the matrix B. Note that Eb' is easily computable 
from Pb'- The algorithm processes columns from 
left to right; for each column, other columns are 
added from the left until a new pivot index appears 
or the column becomes zero. 

Algorithm 2 Reduction Algorithm 

procedure PERSISTENCELEFTRlGHT(i?) 

[0,...,0] 

for j = 1,..., AT do 

while Aj + 0 and L[pivot(Aj)] 0 do 
k <- pivot(Aj ) 

I <- L[pivot(Aj)] 

^3 ^3 ~ • Ai 

if Aj + 0 then 

i t- pivot (Aj) 

L[i\ ^ j 

return {{i,j) \ L[i] = j and j + 0} 


PersistenceData. The basis {wi,... ,wk} rep¬ 
resents the decomposition of the persistent homol¬ 
ogy, and hence, we obtain the persistence data from 
Pb' as follows: 


with 






where {(A, ji), ■ • • , (*s, js)} = {(bj) « Pr I dimcTi = 
I and {rci,rcr^) € C/J, t/^ = {p e | d^{p,A) > e}, 
and {zs+i,-",Zs+t} = {i ^ Er \ dima^ = i}. Hence, 
we have m = 2s + t. The condition e [4 

guarantees the uniqueness of m from Lemma |3.1| 
The following lemma is derived easily from the 
explicit form of h{r). 

Lemma 3.10. The map h is of class C°° on 


We remark that it is sufficient to treat the i-1, 
i, and the {i + l)-simplices in Algorithm if we 
want the persistence diagram for a single value 
of i. 


3.3.4- The map f 

It follows from the chain rule Df{u) = Dh{r) o 
Dg{u) that the explicit form of the derivative 
Df{u) is given by 


with 


Df{u) 


Dfdniu) 

DfiM, ’ 


Dfdniu) = 


dr „ 


du 


du 


du 

<J A 


du 


Dfinfiu) = 


dPa 


du 


dr^ 


du 


), for (i, j) e Pb' and dim cji = £, 
(vcr-, 00 ), for i € Eb' and dima^ = £. 


Note that the persistence data is independent of the 
choice of algorithm from the unique decomposition 
of persistent homology. 

As a result of the above algorithms, the map h : 
is expressed by 


h{r) 


^fin(F) 
hinfic ) 


Proposition 3.11. The map f is of elass on 
OvR OLRd Oaip- Moreover, the derivatives are inde¬ 
pendent of the choiee of algorithms up to permuta¬ 
tions of coordinates in K™. 


Proof. The first statement follows from the chain 
rule and Lemmas 3.6 3.8 and 3.10 For the second 
statement, let us assume two different expressions 
h and h'. From the uniqueness of the persistence 
data, we can express h and h' by appropriate per- 
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mutations if necessary as 
h{u) = 


hfin{u) 

hinf{u) 






and 


h'(u) = 


^Sn(u) = 


Kn(u) 

Kn{(u) 

T<y,, 

1 


j 


. Mnf(w) = 


with 




’"a-,/ 


for all k. On the other hand, it follows from the def¬ 
initions of the Vietoris-Rips and the alpha general 
positions of u that there nniquely exist the attach¬ 
ing simplices rjk and such that 

, ^ <7 / 

for each k. Hence, this leads to 

dn. ^'■"4 


du du du ’ du du du 

and completes the proof of the second statement. 

□ 


Remark 3.12. Since (3.2) and (3.3) are 1-Lipschitz, 


it is reasonable to have differentiability from 
Rademacher’s theorem m- The discussion in this 
section provides a constructive proof of this fact, 
and furthermore, shows that the derivatives are in¬ 
dependent of the choice of algorithm. 


4. Continuation 

4.I. Continuation by Newton-Raphson Method 
We first recall the standard Newton-Raphson 
continuation method m- Let U be an open set in 
M” and (^ : {7 x M -s- K” be a mapping. Suppose 
that (u, A) 6 17 X M satisfies (p{u, A) = 0. Our pur¬ 
pose is to solve ip(u, A) = 0 with respect to ueU for 


a given A. The existence and the local uniqueness 
of the solution u = u\ is guaranteed by the implicit 
function theorem when Du<f{u, A) is regular and A 
is sufficiently close to A. 

In practical computations, we hnd the solution 
u\ for each A by iteratively solving linear equations 
as follows. By taking an appropriate initial point 
the linear approximation of ip{u,X) at each 
iteration step j is given by 

(p{u, A) « (p{u^^\X) + Du(p{u^^\ X){u - 
Setting the right hand side to be zero 

ip{u^^\X) + Du(p{u^^^ ,X){u- = 0, 

we obtain an update of the approximate solution 

ip(tJ'^\X). (4.1) 

This iteration method is called the Newton- 
Raphson Method, and the convergence of the it¬ 
erations under suitable regularity of derivatives is 
well studied (e.g., in])- 

The continuation of the solution (m. A) to a pa¬ 
rameter X' is achieved by gradually changing the 
parameter A. That is, for Aq = A < Ai < ••• < Xn = X', 
the Newton-Raphson method is applied for each pa¬ 
rameter Xi by setting = u\._^ with u\g = u 


4-2. Newton-Raphson Method by Pseudo-Inverse 
Let f - U ^ K"* be a persistence map defined on 
an open set U c K". Let us define a map 


F:Ux] 


(4.2) 


by F(u,v) = f{u) - V. For a given pair (4,7) satis¬ 
fying F(u,v) = 0 and v close to v, our interest is to 
solve F{u,v) = 0 with respect to u eU. The exis¬ 
tence of the solution is guaranteed when Df{u) is 
surjective and v is sufficiently close to v. Hence, for 
the rest of the paper, we add the assumption that 
m < n. In this case, the solutions form an n - m 
dimensional manifold. 

The basic strategy to solve F{u, u) = 0 is the same 


as in Section 4.1 and we derive an iteration method 


for , J = 0,1,..., IV, with = u converging to 
the solution. Namely, the linear approximation of 
F at leads to 


+ F{u^^\ v) = 0, (4.3) 

where = DuF{u^^\v), and we derive an 

update by solving the linear equation with 
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respect to u. However, we note that the linear equa¬ 


tion (4.3) is defined having domain and image with 
different dimensions in general. Thus is 

an m X n rectangular matrix with m < n. 

In general, let us consider a linear equation 


Ax = 6, A e M„ 


be] 


(4.4) 


with m < n. For solving this type of linear equa¬ 
tions, we first recall the concept of pseudo-inverse 
and explain its relation to solutions of the linear 


system (4.4). 


For a matrix A 6 


), there exists a unique 


matrix X e satisfying the so-called Pen¬ 

rose equations: 


AX A = A, XAX = X, 

(AX)'^ = AX, (XAf = XA, 


where is the transpose matrix of A. The unique 
matrix solution X is called the pseudo-inverse of A 
and denoted by Hb 

An explicit formula to construct A^ is given 
for instance in na [n]. Assume that the matrix 
A € has the rank k < m. Then, it has 

a singular value decomposition (SVD) of the form 
A = where V and W are orthogonal ma¬ 

trices, and the matrix E = (crij) e Mm,n(^) has 
CTij- = 0 for all i + j, and > 172,2 ^ ”• ^ o’k,k > 
CTfc+i,fc+i = ••• = (Tm.m = 0. The numbers at := (Ti^i, 
i = are called the singular values of the 

matrix A. From the SVD of the matrix A, the 
pseudo-inverse can be obtained by the formula 


At = ITEtv^, 


where E^ = {(jj j) is the nx m matrix with aj j = 0 
for all * j, (jt. = for f = 1,..., fc, and o-t. = 0 


for k < i < m. 


Equation (4.4) has a solution for x if and only if 
b 6 im(A). In such a case, there is a unique min¬ 
imum norm solution x of (4.4), meaning that the 
Euclidean norm ||x|| is minimum among all the so¬ 
lutions of (4.4). If & ^ im(A), then a least-squares 
solution X of (4.4) is a vector minimizing the error 
||A 2 :-&||. The following proposition provides us with 
relations between the pseudo-inverse and solutions 
of the linear equation. Eor a proof, the reader may 
refer to US] for instance. 


Proposition 4.1. Let Ae andbeW^. If 

b e im(A), then the unique minimum norm solution 
of Ax = b is given by x = A^b. If b i im(A), among 


the least-squares solutions of Ax = b, A^b is the one 
of minimum norm. 


This proposition provides us with a method for 
finding the minimum norm solution of the equation 


(4.3). Namely, we update the approximate solution 
^u/by 


u 


0 +1) _ yO) _ , v ), 


(4.5) 


where Aj = Df{m^'). Note that, from the mini¬ 
mality condition on the norm, the update is 

chosen to be closest to and this is a natural 
choice for the purpose of continuation. 


The convergence of the iterations (4.5) is stud¬ 
ied in [20] (see also [2T|), and we summarize it as 
follows. 


Proposition 4.2. Let r > 0 be a positive real num¬ 
ber such that f 6 C^(Br(u)). Let a, (3 be posi¬ 
tive constants such that, for all u,w e Br{u) with 
u - w e iml?/(iy)^, the followings hold: 


\\Df{w){u -w)- f{u) + f{w)\\ < a\\u - w\ 
\\{Df{w) ^ - Df{u) 0/(m)|| < (3\\u - wll, 
Q!||D/( 2 ;)^|| + ,5 = 7 < I ioT all z e Br(u), 
\\Df{uy\\-\\f{u)\\ < 


Then the iteration (4.5) converges to a solution of 


Df{u)'^F{u,v) = 0 


which lies in Bj.{u). 


When m = n, this proposition provides a criterion 
for the convergence of the Newton-Raphson method 

For rank(ZI/(M)) = m, the iteration converges 
to a solution of F{u,v) = 0. On the other hand, 
when rank(T)/(u)) < m, the convergent point u 
does not necessarily satisfy F{u,v) = 0, but only 
implies F{u,v) 6 herDf(u). Thus, in our contin¬ 
uation method, we suppose that all the singular 
values (Ti,..., ajn are positive. 

4 . 3 . Continuation of Point Clouds 


We use the iteration (4.5) for continuations of a 
point cloud. Let f :U ^ be a persistence map. 
Suppose that {us,Vs) e U x K"* is a pair satisfying 
f{us) = Vs- This pair can be regarded as the initial 
point of the continuation. Our task is to continuate 
it to a target persistence data Vt € and obtain 
ut e M" satisfying f{ut)=Vt. 
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As the simplest way of the continuation, we di¬ 
vide the line between Vs and Vt equally into N small 
segments, and for each v = Vg + kAv, k = 1,..., A^, 
where Av = , we apply the iteration method 


N 


(4.5) and obtain the point cloud u satisfying f{u) = 
V. This process is summarized in Algorithm]^ 


Algorithm 3 Continuation of a point cloud 
input Us,Vs,Vt,N 
Av = (vt - Vs)/N, u = Us 
for k = \ N Ao 
V = Vg + kAv 

if PseudoInverseNewton(u,'!;) converges 
u' =PseudoInverseNewton(u, u) 
u = u' 
else 
break 
endif 

output Ut = u 


We can, of course, adaptively choose the length 
of Av at each continuation step. Furthermore, we 
can also adopt an appropriate curve connecting Vg 
to Vt if necessary. 

We also note that the image of / is not generally 
equal to the target space K"*. Hence, if we choose 
Vt in K"* \ im(/), the continuation fails at some 


V = Vg + kAv. See Section 5.2 for such an example. 


It is often the case in practical problems that 
the point clouds need to satisfy some constraints 
gi{u) = 0, i = l,...,r (e.g., conservation laws in 
mechanics). In such a case, we need to solve the 
continuation under these constraints, and its modi¬ 
fication is straightforward. Namely, we first extend 


the original setting (4.2) to 


F : C7 X M”" 

by F{u,v) = {f{u)-v,gj(jA,...,gr{u)). Then, by 


replacing Aj and F in (4.5) with Aj = D^F and 


F, respectively, we obtain the appropriate formulas 
for continuation of point clouds under constraints 
gi{u) = 0, i= l,...,r. 


4.4- Symmetry 

It should be noted that we need to remove sym¬ 
metries induced by translations and rotations in or¬ 
der to isolate a solution u of f{u) = v. For a point 
cloud P in the following restrictions will remove 
these symmetries: 

(i) Fix one of the vertices in P at the origin of 
say uo = (0,0,0). 


(ii) Select one of the vertices in P \ {uo}j say ui. 
This vertex is supposed to move only on the 
line defined by u^ui during the continuation. 

(hi) Select one of the vertices in P \ {uo,ui}, say 
U 2 . This vertex is supposed to move only on 
the plane dehned by uqUi,uqU 2 during the con¬ 
tinuation. 

The first restriction eliminates translation symme¬ 
tries and the second and third restriction eliminate 
rotation symmetries. 

In practice we choose our basis of the coordinate 
system (x, y, z) € in such a way that, in addition 
to Uq being fixed at the origin, we have that Ui 
stays on the x-axis and U 2 stays on the (x, 2 /)-plane 
during the continuation. Hence, in the following, 
let us redefine n = 3M - 6 and 

U= (ui,U2,U3,...,Um), (4.6) 

where ui = xi,U 2 = (x 2 ,y 2 ),Ui = {xi,y^,Zi) for 3 < 
i < M. 

For a general point cloud in K^, these restric¬ 
tion to eliminate the symmetries should be appro¬ 
priately modified. 


4-5. Non-convergence and Zero Singular Values 
The Newton-Raphson method by pseudo-inverse 
does not work well when the matrix Aj in (4.5) has 
a singular value close to zero. In this section, we 
discuss those cases in the alpha and Vietoris-Rips 
filtrations. 


4-5.1. For the case of alpha filtrations 

For the case of alpha filtrations, we show that a 
singular value close to zero appears when there is a 
birth-death pair (6, d) with 6 rs d in the persistence 
diagram P£(Alp(P)). Here, we impose the natural 
assumptions that a point cloud P in [L = 2,3) 
satisfies the general position assumption and con¬ 
sists of at least L + 1 points. 

Let Alp(P, a+) be the alpha complex constructed 
by simplices whose birth radii are smaller than 
or equal to a. Notice that this is different from 
Alp(P, a), where a is equal to the birth radius of a 
simplex. 

Proposition 4.3. Let P be a point cloud and let 
{b,d) be a birth-death pair in Di(Abp{P)) (i=l,2). 
If the birth radius of any simplex is not contained 
in the open interval {b,d), there exist an attaching 
(.-simplex a and an attaching {(+ 1)-simplex r such 
that r^ = b, Ct- = d, and t a. 
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To prove the proposition, we need the following 
lemma. 

Lemma 4.4. Let a be a simplex and r be its at¬ 
taching simplex with a + t. Then, the inclusion 
from |Alp(P,rT-)! to |Alp(P, rT-+)| is a deformation 
retract. In particular, r,- is neither a birth time nor 
a death time. 

To prove the lemma, we recall some properties 
about Voronoi decomposition [22] • For fc + 1 points 
{uq, ... ,Uk} in P, let V {uq, •••, Uk) be defined as 

V{uo,---,Uk) 

={x €M.^ \ d{x,uo) = ■■■ = d(x,Uk), and 

d{x, Uq) < d(x, u) for Vue F\{uo ,.. •, itfe}}. 

For fc = 0, the set V(uo) is called a Voronoi cell. 
From the theory of Voronoi decompositions and De¬ 
launay triangulations, we have the following facts: 

1. V(uo,---,Uk) = V(uo)n---nV(uk). 

2. {uo,---,Uk} forms a Delaunay fc-simplex if and 
only if V(mo, Mfc) is not empty. 

3. V{uo,---,Uk) is closed, convex, and con¬ 
tained in a {L - fc)-dimensional affine sub¬ 
space W{uo,---,Uk) = {a; 6 I d{x,uo) = 

■■■ = d{x,Uk)}. W{uo,---,Uk) is orthogonal to 
the fc-dimensional affine subspace spanned by 

{ziQ 5 ■ ■ • 5 } ■ 

4. The boundary of V{uo,---,Uk) relative to 
W{uQ,---,Uk) is [JviY{uo,-,u^)y{uo,---,Uk,v), 
where Y{uo,---,Uk) = {v e P\{uo,---,Uk} \ 
V{uo,—,Uk) nV{v) + 0}. 

5. If V{uo,---,Uk) is not empty, Y{uo,---,Uk) is 
not empty. 


E = V{ui,U 2 ,uo) from the fact Then, from the 
definition of P = V{ui,U 2 ,uq) and the general po¬ 
sition assumption, the following holds: 

d{E,ui) = d{E,U 2 ) = d{E,uo), 
d{E, u) > d(E, Ml) for Vu € P\{mo, mi, M 2 }. 

We take a new orthogonal coordinate system on 
satisfying 

• P= (0,0) 

• V(mi,M 2 ) is contained in{(t, 0 )|t< 0 }. 

In this coordinate, mq, mi, and M 2 , are described as 

Mo = (rcos0o,rsin6*o), 

Ml = (r cos 01 , r sin 6 * 1 ), 

M 2 = (r COS 01 ,-r sin0i), 

since mo,mi,M2 lie on the same circle whose cen¬ 
ter is E and a is orthogonal to V(mi,M 2 ) (Fig¬ 
ure [^. It follows from a nV{ui,U 2 ) = 0 that 0i 
must be contained in [ 0 , 7 r/ 2 ). Since = (-e,0) € 
V(mi,M2)\V(mi,M2,mo) for small e > 0 , we have 
d(Pe,Mi) < fi(Pe,Mo) and 

0< (i(Pe, Mo)^ - d(Pe, Mi)^ 

= ((r cos 00 + e)^ + (?’sin0o)^) 

- ((rcos 0 i + e)^ + (rsin 0 i)^) 

= 2 re(cos 0 o - cos 0 i). 

Hence - 7 r /2 < -0i < 0o < 0i < 7 r /2 holds and r is an 
obtuse triangle whose longest edge is a. Therefore, 
for a' = mqMi and a" = U 0 U 2 , we can show that 
Tct' < Tt and To-" < Xr from the general position 
assumption and the definition of the birth radius of 
a simplex oj 


6. A Delaunay fc-simplex {Mo,...,Mfc} is attach¬ 
ing if and only if this simplex and V (mo, •••, Uk) 
have a non-empty intersection. 


Proof. (Lemma 4.4) First we consider the case of 
L = 2,dimer = 1, and dimr = 2. Let mi,M 2 be 
the endpoints of a and mq be the other vertex of 
T. From the facts and above, V(mi,M2) is 
not empty and contained in the perpendicular bi¬ 
sector of cr. From the fact V(mi,M 2 ) is not 
empty and V(mi,M 2 ) is a line with one endpoint 
or with two endpoints. Since a is not attaching, 
V(mi,M 2 ) n ct = 0 from the fact Let P be the 
endpoint of V(mi,M 2 ) close to cr, which is given by 


ruj = infjthe radius of P | P is an empty ball 
satisfying dB d (vertices of w)}. 

Since any other birth radii are different from 
Xr by the general position assumption, we have 
Alp(P, r,- + )\Alp(P, r,-) = {cr, r}. Therefore, the in¬ 
clusion 

|{cr',cr",Mo,Mi,M2}| ^ |{r, cr, cr', cr", Mq, Mi, M2, }|, 

leads to the desired deformation retract. 

The case of P = 3,dimer = 1, and dimr = 2 can 
be proved in the same way by considering the plane 
that contains three vertices of r. 
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Figure 2: Orthogonal coordinate system for an attaching 2- 
simplex r = {uq i “1 > ^^2 } • 


The case of L = 3,dimtT = 2, and dimr = 3 
can also be proved in a similar way by replacing 
V(ui,U 2 ) in the above argument to V{ui,U 2 tU^)- 
In this case, the 3-simplex r is given by a tetrahe¬ 
dron whose four vertices are on a half side of the 
circumsphere of the tetrahedron. □ 

The following corollary is a straightforward con¬ 
sequence of Lemma |4.4| 

Corollary 4.5. For a simplex t, if is either a 
birth time or a death time in the persistence dia¬ 
gram, T is an attaching simplex and does not attach 
any faces ofr. 

Proof. (Proposition |4.3[ ) Since d is a death time in 
Z)^(Alp(P)), an (£+ l)-simplex r is born at time d, 
hence d = r^.. From Corollary |4.5[ r is an attaching 
simplex and does not attach any faces of t. With 
the general position assumption, t is the unique 
simplex satisfying d = r,-. 

Next, we consider the change at time b. Since 
the birth radius of any simplex is not contained 
in (6, d) and r is unique as above, Alp(P, a) = 
Alp(P,rT-+)\{T} for any a e (6,d]. Therefore, one 
of the ^-dimensional faces tr of r appears at time 
b. Otherwise, the generator by dr keeps unchanged 
at time b and this contradicts the fact that b is the 
birth time of the pair {b,d). The simplex cr is an at¬ 
taching simplex from Corollary |4.5[ and hence this 
concludes the proof. □ 

Now, we show that if a birth-death pair (&, d) is 
close enough b d for the birth radius of any sim¬ 
plex not to be contained in {b,d), then a singular 
value close to zero appears in the derivative of the 
persistence map. In this case, from Proposition 
there exist an (€+1) -simplex t and its face a satisfy¬ 
ing r^ = d and r^ = b. Let pi,... ,pi be the vertices 
of a and po,pi,... ,pi be the vertices of r. Let cq 
be the center of the minimal circumsphere of t and 
Cl be the center of the minimal circumsphere of cr. 


Let r be the length of PqCq, the radius of minimal 
circumsphere of t, and 0 be the angle of ^copici 
(Figure]^. From the factof the Voronoi decom¬ 
position, cgCi is orthogonal to tr. Therefore we can 
represent r,- and ra- as follows: 

Cr = c, 

To- = rcosd. 




Figure 3: In the left figure, cq is the center of the circum- 
circle of the triangle {p 0 iPiiP 2 }, ci is the center of pip 2 i t' 
is the radius of the circumcircle, and 6 is ^copici. In the 
right figure, cq is the center of the circumsphere of the tetra¬ 
hedron {po,pi,p 2 ,P 3 }, Cl it the center of the circumcircle of 
the triangle {pi,P 2 ^P 3 }^ ^ is the radius of the circumsphere, 
and 6 is ^copici. 


We consider a map f-U 


given by 


f :u< 


assigning one birth-death pair. This map is one 
component of the persistence map and can be de¬ 
composed as follows 


u h 




(r,6») h 


h 


{rr,ra) ■ 


We can easily show the following facts: 


• r^ = r„ if and only if 0 = 0. 


• If 0 = 0, 


holds. 


Df2 


1 0 
1 0 


From these facts, we have 0 « 0 for r,- ^ r^ and 


Df2 


1 0 
1 0 


Since the matrix of the right hand side is not sur¬ 
jective, Df = Df 2 ° Dfi has a singular value close 
to zero, and hence the derivative of the persistence 
map has a singular value close to zero. 
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4-5.2. For the case of Vietoris-Rips filtrations 
For the case of Vietoris-Rips filtrations, a birth- 
death pair {b,d) with b d does not necessarily 
imply the existence of a singular value close to zero. 
We show such an example for a point cloud P = 
{A, B,C, D} in and Di(VR{P)). We assume 
the followings: 

• |AB|, |RC|, \CD\ < \AD\ < \BD\ < \AC\. 

• |AD| Ri \BD\. 

From the assumption, the triangle ABD is close 
to an isosceles triangle. The persistence diagram 
Di(VR(P)) has a unique birth-death pair (&,d) 
where b = \AD\ and d = \BD\ (Figure]^. 



Figure 4: A point cloud with four points. Di(VR(P)) has 
a unique birth-death pair (6, d) 


The persistence map / is 



and from the computation in Lemma |3.6[ Df • 
(Df)^ is 




2 cosd 
cos 0 2 ’ 


where 0 = ^ ADB. Since eigenvalues of Df ■ (Df)"^ 
are squares of singular values of Df, the singular 
values are \/2 ± cos 0 and are away from zero, al¬ 
though 6 Ri d. 


5. Computations 

In this section, we present some numerical exam¬ 
ples of continuations of point clouds using persis¬ 
tence diagrams. Alpha filtrations are used for all 
examples, and the coordinate system described in 
Section is adopted to eliminate the translation 
and rotation symmetries. 

Example 5.1 (Deformation of a tetrahedron). As 
a first example, we consider alpha filtrations con¬ 
structed from four points in (a tetrahedron) 
and apply our continuation algorithm to D 2 . We 


take as the initial point cloud P = {uo,ui,U 2 ,U 3 }, 
where uq = (0,0,0), ui = (8,0,0), U 2 = (5,6,0), 
and U 3 = (4,2,6). The 2nd persistence diagram 
of P is D 2 = {(4.42719,4.59015)}. From this ini¬ 
tial data, we try to deform D 2 to the target persis¬ 
tence data {(8.42719,8.89015)} using our continua¬ 
tion method. Due to the coordinate system adopted 


to eliminate symmetries, (4.6), we have that the de¬ 


gree of freedom of the point cloud is six, and the 
persistence map can be expressed as /:t/ c K® ^ 
For this example we used ||Au|| = 0.01 as the 
step size in the continuation, and e = Cjfj 

Figures and show the point clouds and 
the persistence diagrams during the continua¬ 
tion process. In this computation we suc¬ 
cessfully reached the target persistence diagram 
{(8.42719,8.89015)}. However notice that there 
seems to be a non-smooth point on the blue curves 
given by the point clouds during continuation in 
Figure To try to understand this event, we 
look at the birth radii of the 2-simplices during 
the continuation process (Figure [^. Notice that at 
some point during the continuation process, two of 
the birth radii coincide, hence breaking the second 
condition of the alpha general position assumption 
(Definition |2.2[ ). This point corresponds exactly to 
the non-smooth point in Figure This is a point 
where the derivative Df is not uniquely defined, 
hence the continuation curve is not smooth there. 
The fact that we are performing the continuation 
numerically, and hence the birth radii are very close 
but not exactly equal, makes it possible to go for¬ 
ward with the continuation process. At each step 
of the continuation the birth radius that is slightly 
larger is used to compute Df, and the 2-simplex 
corresponding to this radius can change at each con¬ 
tinuation step. We can see this from the singular 
values of Df show in Figure 

Example 5.2 (Image of the persistence map). In 
this example we use the same point cloud (tetra¬ 


hedron) as in Example 5.1 to explore the image of 
the persistence map. For a tetrahedron, the image 
of the persistence map is the strip region between 
the diagonal and the line of persistence diagrams of 
regular tetrahedrons (see Figurej^and Theorem 6.1 
in the Appendix for a proof of the Di case). In this 
example we used || Au|| = 0.001 and e = 0. 


^The 2-dimensional diagram D 2 of a tetrahedron has at 
most one birth-death pair, hence we cannot cut off points 
close to the diagonal and therefore we use e = 0. 
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Figure 5: Deformation of the point cloud during continua¬ 
tion. The black dots are the initial point cloud and the red 
dots are the final one. The blue curves show the movement 
of the points during the continuation. Note that the first ver¬ 
tex is fixed at the origin and the second one is only allowed 
to move along the x-axis. 



Figure 7: Birth radii of the 2-simplices along the continua¬ 
tion curves in Figure 
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Figure 6: Persistence diagrams during the continuation 

process. The black dot represents the initial persistence di- Figure 8: Singular values of the derivative Df along the 

agram and the red dot represents the target diagram. The continuation in Figure]^ 

blue dots represent the persistence diagrams during the con¬ 
tinuation process. 
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From the initial persistence diagram D 2 = 
{(4.42719,4.59015)} we try to continue to the tar¬ 
get persistence {(6.42719,7.09015)}, which is out¬ 
side of the image of the persistence map. Hence, as 
expected, in this case we fail to reach the target per¬ 
sistence and can only continue up to the boundary 
of the image, as we can see in Figure As we ap¬ 
proach the boundary of the image, the method fails 
because the number of Newton iterations needed 
for convergence increases dramatically as is shown 
in Figure In the last steps of the continuation 
the birth radii of the 1 -simplices are very similar 
and the birth radii of the 2 -simplices are all virtu¬ 
ally the same as we can see in Figures [m and [T^ 
hence confirming that we have continued to a regu¬ 
lar tetrahedron. Notice also from Figures [TT] and [12] 
that two or more birth radii are equal during the 
continuation, hence the general position assump¬ 
tion (Definition |2.2[ ) is violated. However, as noted 
in Example 1 5.1 1 the continuation method still works 
as long as we are within the image of the persistence 
map. 



Figure 9: The image of the persistence map of D 2 for one 
tetrahedron (four points in R^) is the shaded region between 
the diagonal and the line Li of the persistence diagrams of 
regular tetrahedrons (obtained by similarity deformations). 
The lines L 2 and L 3 correspond to the persistence diagrams 
of tetrahedrons with three and two congruent faces, respec¬ 
tively. The blue line shows the continuation of the initial 
point cloud all the way to the boundary of the image. 

Example 5.3 (Towards the diagonal). 
In this example we take the point cloud 
P = {uo,ui,U 2 ,U 3 }, where uq = (0,0,0), 

Ml = (9.991,0,0), m 2 = (4.9955,8.65246,0), and 
M 3 = (4.9955,2.88415,8.15762), which represents 



Figure 10: Number Newton iterations during each step of 
the continuation in Example] 5.2 1 



Figure 11: Birth radii of the 1-simplices along the continua¬ 
tion in Example 1 5.2 1 



Figure 12: Birth radii of the 2-simplices along the continua¬ 
tion in Example 1 5.2 1 
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a nearly regular tetrahedron, and try to continue 
the persistence diagram D 2 towards the diagonal. 
From the discussion in Section 14.51 this will lead 
to singular values close to zero. The 2nd persis¬ 
tence diagram of P is D 2 = {(5.76831,6.11821)}, 
and the target persistence diagram is set to be 
{(5.94841,5.94841)} which is on the diagonal. In 
the computations we used || Au|| = 0.001 and e = 0. 

In this case the computations work well all the 
way to a point nearly on the diagonal. In Figure [13] 
we show the persistence diagrams along the contin¬ 
uation. The singular values of the derivative are 
shown in Figure 14 As expected, one of the singu¬ 
lar values approaches zero towards the end of the 
continuation. However, in spit of this, the contin¬ 
uation works all the way to a point essentially on 
the diagonal. Note that we cannot continue to a 
point exactly on the diagonal, since the persistence 
diagram would be empty in that case. However the 
persistence diagram that we arrive at the end of the 
continuation in this example is only “numerically” 
on the diagonal, that is, it is on the diagonal up to 
the error tolerance of the Newton-Raphson method. 
As described in Section 14.21 the method will fail if 
we try to continue to a point exactly on the diago¬ 
nal, since in that case we would have a zero singular 
value. 



Figure 13: Persistence diagrams D 2 during the continuation 
process. The black dot represents the initial persistence di¬ 
agram and the red dot represents the final diagram. Notice 
that the final diagram is “numerically” on the diagonal. 



Figure 14: Singular values of the derivative Df along the 
continuation in Example |5.3| Notice that one of the singular 
values approach zero as we reach the end of the continuation. 

Di = {(0.758288,0.803195), (0.776209,0.834393)} 
to the target 1 -dimensional persistence diagram 
{(0.770801,0.817236), (0.798346,0.863075)}. In 
these computations we used || Aujl = 0.001 and e = 0 . 
There were only two points in the persistence dia¬ 
gram during all the steps of the continuation, hence 
the choice e = 0. The computations worked well all 
the way to the target persistence. In Figure [Ts] we 
show the point cloud and the persistence diagrams 
along the continuation. 



Figure 15: Persistence diagrams Di during the continuation 
process. The black dots represent the initial persistence di¬ 
agram and the red dots represent the final diagram. 


Example 5.4 (Continuation of Di). In this 
example we take the point cloud data P = 
{uo,ui,U 2 ,U 3 }, where uq = ( 0 , 0 , 0 ), ui = ( 1 , 0 , 0 ), 
U 2 = (1.1,1.2,0), and M 3 = (0.5,0.6,1.3), and 
try to continue the initial persistence diagram 


Example 5.5 (Deformation of a dodecahedron). 
In this example we take as the initial point cloud 
the vertices of a regular dodecahedron and try to 
apply our continuation method to increase both the 
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birth and the death radii of the D 2 generator. In 
these computations we used ||Au|| = 0.01 and e = 
10“^. The continuation works all the way to the 
target persistence diagram. Figure shows the 
deformation of the point cloud and Figure [T7| shows 
the diagrams along the continuation. 



Figure 16: Deformation of the point cloud during continua¬ 
tion. The black dots are the initial point cloud and the red 
dots are the final one. The blue curves show the movement 
of the points during the continuation. 


Example 5.6 (Deformation of a sphere). In this 
example we take as the initial point cloud 100 uni¬ 
form point on a sphere and try to apply our con¬ 
tinuation method to the largest generator of II 2 . 
In these computations we used ||Au|| = 0.03 and 
e = 10~®. The continuation works all the way to the 
target persistence diagram. Figureshows the de¬ 
formation of the point cloud and Figure [T^ shows 
the diagrams along the continuation. 



Figure 17: Persistence diagrams D 2 during the continuation 
process. The black dot represents the initial persistence di¬ 
agram and the red dot represents the final diagram. 



Figure 18: Deformation of the point cloud during continua¬ 
tion. The black dots are the initial point cloud and the red 
dots are the final one. The blue dots show the movement of 
the points during the continuation. 
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Birth 

Figure 19: Persistence diagrams D 2 during the continuation 
process. The black dot represents the initial persistence di¬ 
agram and the red dot represents the final diagram. 

6. Conclusion 

In this paper, we have studied the continuation 
of point clouds by persistence diagrams. In the 
following, we list some future improvement of our 
method. 

1. In the presented method, we have treated per¬ 
sistence diagrams in an Euclidean space whose 
dimension is determined by the input persis¬ 
tence diagram. This vectorization is simple 
and describes the essential part for the con¬ 
tinuation method. However, it does not allow 
to change the number of generators in the per¬ 
sistence diagrams during the continuation be¬ 
cause of the fixed dimension of the Euclidean 
space. To overcome this restriction, it can 
be useful to vectorize the persistence diagrams 
into a bigger space and construct a similar con¬ 
tinuation method. The space of persistence 
landscape |24j or a vectorization using kernel 
methods [5S] should be considered as possible 
candidates. 

2. Our algorithm for computing persistence dia¬ 
grams in this paper is not sophisticated, and 
hence, there is room for improvement. Stan¬ 
dard reduction methods such as [55] can be im¬ 
plemented and will reduce the computational 
cost. Furthermore, since the changes in the 
point clouds at each step of the continuation is 
supposed to be small, the vineyard algorithm 
m can effectively work for fast computations. 
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Appendix. Image of the Persistence Map of 
ZDi for a triangle 

Theorem 6.1. If a point cloud P has only three 
points, Di(Alp{P)) has at most one hirth-death 
pair. //ZIi(Alp(P)) has such a pair {b,d), the ra¬ 
tio d/b is smaller than or equal to 2/\/3. Moreover, 
djb = 2/\/3 if and only if the triangle is regular. 

Proof. From the basic properties of alpha filtra- 
tions, IIi(Alp(P)) = {{b,d)} if and only if the tri¬ 
angle is acute. If not, £)i(Alp(P)) is empty. Hence 
we assume that the triangle is acute. 

Let Po,Pi,P 2 be the three vertices of the triangle, 
c be the center of the circumcircle, r be the radius of 
the circumcircle, and 0o,Oi,92 be ^cpip 2 , ^cp 2 Po, 
and ^cpopi, respectively (Figure]^. 



Figure 20: The relation of the circumcircle and internal an¬ 
gles of the triangle {po,Pi,P 2 } 

Hence, \piP 2 \ = 2rcos6»o, \p 2 Po\ = 2rcos6»i, \popi \ = 
2r cos 02 and the birth time b is 

b = max{rcos0o,^cos0i,rcos02}, 

and the death time d is r. The ratio of d/b is 

r/maxjrcos0o,rcos 0i,rcos 02 } 

= (max{ cos 0o: cos , cos 02 }) ~ ^ ■ 

Hence the problem is minimizing 

max{ cos 0o j cos 0i, cos 02 } 
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subject to 

> 0 (the triangle is acute) and 
9o + 01 + 02 = 7r/2 (sum of the internal angles). 

Since cos 0 is monotonely decreasing on the interval 

[0,7r/2], a minimum attains at 0q = 9i = 02 = tt/G 

and the minimum is cos(7r/6) = \/3/2. □ 
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